#include<../data_tools/usefull_routines.c>

int main(int argc, char *argv[])
{
  int i,j,npart,nbines,ndisk,limpart,counter;
  double *radio,rmin,rmax,deltar,rini,rfin,rpot;
  double *potential,meanpot,pot[6];
  size_t *ind;
  char *infile,outfiles[1000];
  particulas *diskpart;
  infile = argv[1];
  limpart = 100;

  FILE *diskfile;
  diskfile = fopen(infile,"r");
  ndisk = contadorParticulas(infile);

  diskpart = (particulas *)malloc((size_t)ndisk*sizeof(particulas));
  if(diskpart == NULL){
    printf("Allocation error routines:81\n");
    exit(0);
  }

  radio = (double *)malloc((size_t)ndisk*sizeof(double));
  if(radio == NULL){
    printf("Allocation error routines:81\n");
    exit(0);
  }

  potential = (double *)malloc((size_t)ndisk*sizeof(double));
  if(radio == NULL){
    printf("Allocation error routines:81\n");
    exit(0);
  }

  ind = (size_t *)malloc((size_t)ndisk*sizeof(size_t));
  if(ind == NULL){
    printf("Allocation error routines:81\n");
    exit(0);
  }

  for(i=0; i<ndisk; i++)
    {
      fscanf(diskfile,"%d %e %e %e %e %e %e %e %e",
	     &diskpart[i].id,
	     &diskpart[i].pos[0],&diskpart[i].pos[1],&diskpart[i].pos[2],
	     &diskpart[i].vel[0],&diskpart[i].vel[1],&diskpart[i].vel[2],
	     &diskpart[i].mass,
	     &diskpart[i].pot);
      /*      radio[i] = sqrt(diskpart[i].pos[0]*diskpart[i].pos[0]+
		      diskpart[i].pos[1]*diskpart[i].pos[1]+
		      diskpart[i].pos[2]*diskpart[i].pos[2]);*/
      radio[i] = sqrt(diskpart[i].pos[0]*diskpart[i].pos[0]+
		      diskpart[i].pos[2]*diskpart[i].pos[2]);
      potential[i] = diskpart[i].pot;
      //ind[i] = i;
    }
  
  gsl_sort_index(ind,radio,1,ndisk);
      
  
  rpot = 0.5;
  counter = 0;
  for(i=0; i<ndisk; i++)
    {
      //      if((radio[ind[i]]>=(rpot-0.1)) && (radio[ind[i]]<=(rpot+0.1)))
      if((diskpart[i].pos[0]>=(rpot-0.1)) && (diskpart[i].pos[0]<=(rpot+0.1)))
	{
	  if((diskpart[i].pos[2]>=(0.0-0.1)) && (diskpart[i].pos[2]<=(0.0+0.1)))
	    {
	      pot[counter] = potential[i];
	      //	      printf("%d %e \t %e\n",counter,radio[ind[i]],pot[counter]);
	      printf("%d %e \t %e\n",counter,
		     sqrt(diskpart[i].pos[0]*diskpart[i].pos[0]+diskpart[i].pos[2]*diskpart[i].pos[2]),
		     pot[counter]);
	      rpot = rpot + 0.5; 
	      counter++;
	      if(rpot==3.5) break;
	    }
	}
    }
  /*
  FILE *pf;  
  for(j=0; j<=5; j++)
    {
      sprintf(outfiles,"equal_pot_%s_%d.dat",infile,j+1);
      //printf("%s\n",outfiles);
      pf = fopen(outfiles,"w");

      for(i=0; i<ndisk; i++)
	{
	  if((potential[i]>=(pot[j]-100.0)) && (potential[i]<=(pot[j]+100.0)))
	    {
	      fprintf(pf,"%e %e %e %e %e\n",
		      diskpart[i].pos[0],
		      diskpart[i].pos[1],
		      diskpart[i].pos[2],
		      pot[j],
		      diskpart[i].pot);     
	    }
	}
      fclose(pf);
    }
*/
  free(diskpart);
  free(radio);
  free(potential);
  free(ind);

  return 0;
}
